"""
    DiagnosticReaders

Read and validate GSI diagnostic files.

This module provides functions to read GSI diagnostic files (both Fortran and Julia
generated) for validation and comparison purposes.
"""
module DiagnosticReaders

using Printf
using Statistics

# Import required modules
include("FortranBinaryIO.jl")
include("DiagnosticFileFormat.jl")

using .FortranBinaryIO
using .DiagnosticFileFormat

export read_diag_obs_file, read_fit_statistics, validate_diagnostic_file
export compare_diagnostic_files, DiagnosticFileInfo

"""
    DiagnosticFileInfo

Information about a diagnostic file.
"""
struct DiagnosticFileInfo
    filename::String
    file_size::Int
    n_records::Int
    n_observations::Int
    var_name::String
    var_type::String
    grid_size::Tuple{Int, Int, Int}
    observations::Vector{ObservationRecord}
end

"""
    read_diag_obs_file(filename::String)

Read a GSI observation diagnostic file (analysis or background).

Returns a DiagnosticFileInfo struct with all file contents.
"""
function read_diag_obs_file(filename::String)

    println("Reading diagnostic file: $filename")

    if !isfile(filename)
        error("File not found: $filename")
    end

    file_size = filesize(filename)
    println("  File size: $file_size bytes")

    ff = open_fortran_file(filename, "r", byte_order=:BigEndian)

    observations = ObservationRecord[]
    var_name = ""
    var_type = ""
    grid_size = (0, 0, 0)
    n_records = 0

    try
        # Read header
        header = read_fortran_header(ff)
        if header !== nothing
            var_name, var_type, metadata = header
            println("  Variable: $var_name ($var_type)")
            n_records += 1
        end

        # Read grid metadata
        grid_data = read_fortran_record(ff, UInt8)
        if grid_data !== nothing
            grid_size = parse_grid_metadata(grid_data)
            println("  Grid size: $grid_size")
            n_records += 1
        end

        # Read observation records
        obs_count = 0
        while true
            obs_data = read_fortran_record(ff, UInt8)
            if obs_data === nothing
                break
            end

            n_records += 1

            # Try to parse as observation records
            # Each observation is approximately 164 bytes
            record_size = length(obs_data)
            obs_record_size = 164  # Approximate size from pack_observation_record

            if record_size >= obs_record_size
                # Parse multiple observations from this record
                n_obs_in_record = record_size ÷ obs_record_size

                for i in 1:n_obs_in_record
                    offset = (i-1) * obs_record_size + 1
                    if offset + obs_record_size - 1 <= record_size
                        try
                            obs_bytes = obs_data[offset:offset+obs_record_size-1]
                            obs = unpack_observation_record(obs_bytes)
                            push!(observations, obs)
                            obs_count += 1
                        catch e
                            # Skip malformed observation
                            @debug "Failed to parse observation $obs_count: $e"
                        end
                    end
                end
            end
        end

        println("  Total records: $n_records")
        println("  Observations read: $obs_count")

    catch e
        println("  Warning: Error reading file: $e")
    finally
        close_fortran_file(ff)
    end

    return DiagnosticFileInfo(
        filename,
        file_size,
        n_records,
        length(observations),
        var_name,
        var_type,
        grid_size,
        observations
    )
end

"""
    parse_grid_metadata(data::Vector{UInt8})

Parse grid metadata from binary data.
"""
function parse_grid_metadata(data::Vector{UInt8})
    if length(data) < 32
        return (0, 0, 0)
    end

    io = IOBuffer(data)

    try
        grid_type = ntoh(read(io, Int32))
        nx = ntoh(read(io, Int32))
        ny = ntoh(read(io, Int32))
        nz = ntoh(read(io, Int32))

        return (Int(nx), Int(ny), Int(nz))
    catch
        return (0, 0, 0)
    end
end

"""
    read_fit_statistics(filename::String)

Read a GSI fit statistics file (fort.*).

Returns a dictionary with statistics.
"""
function read_fit_statistics(filename::String)

    println("Reading fit statistics: $filename")

    if !isfile(filename)
        @warn "File not found: $filename"
        return Dict{String, Float64}()
    end

    stats = Dict{String, Float64}()

    open(filename, "r") do f
        lines = readlines(f)

        for line in lines
            # Parse o-g line
            if occursin("o-g", line) && occursin("asm all", line)
                parts = split(line)
                if length(parts) >= 9
                    try
                        stats["count"] = parse(Float64, parts[6])
                        stats["bias"] = parse(Float64, parts[7])
                        stats["rms"] = parse(Float64, parts[8])
                        stats["cpen"] = parse(Float64, parts[9])
                        if length(parts) >= 10
                            stats["qcpen"] = parse(Float64, parts[10])
                        end
                    catch e
                        @debug "Failed to parse fit statistics line: $e"
                    end
                end
            end

            # Parse rejection count
            if occursin("failed gross test", line)
                m = match(r"=\s*(\d+)", line)
                if m !== nothing
                    stats["n_rejected"] = parse(Float64, m.captures[1])
                end
            end
        end
    end

    println("  Statistics: $stats")

    return stats
end

"""
    validate_diagnostic_file(filename::String)

Validate a diagnostic file structure and contents.

Returns a tuple of (is_valid, messages) where messages contains
validation results.
"""
function validate_diagnostic_file(filename::String)

    messages = String[]
    is_valid = true

    try
        info = read_diag_obs_file(filename)

        # Check file size
        if info.file_size == 0
            push!(messages, "ERROR: File is empty")
            is_valid = false
        else
            push!(messages, "✓ File size: $(info.file_size) bytes")
        end

        # Check records
        if info.n_records == 0
            push!(messages, "ERROR: No records found")
            is_valid = false
        else
            push!(messages, "✓ Records: $(info.n_records)")
        end

        # Check observations
        if info.n_observations == 0
            push!(messages, "WARNING: No observations found")
        else
            push!(messages, "✓ Observations: $(info.n_observations)")
        end

        # Check header
        if isempty(info.var_name)
            push!(messages, "WARNING: No variable name in header")
        else
            push!(messages, "✓ Variable: $(info.var_name)")
        end

        # Check grid size
        if all(info.grid_size .== 0)
            push!(messages, "WARNING: Grid size not found")
        else
            push!(messages, "✓ Grid: $(info.grid_size)")
        end

        # Validate observation values
        if info.n_observations > 0
            valid_lats = count(obs -> -90 <= obs.lat <= 90, info.observations)
            valid_lons = count(obs -> -180 <= obs.lon <= 360, info.observations)

            if valid_lats < info.n_observations
                push!(messages, "WARNING: $(info.n_observations - valid_lats) invalid latitudes")
            else
                push!(messages, "✓ All latitudes valid")
            end

            if valid_lons < info.n_observations
                push!(messages, "WARNING: $(info.n_observations - valid_lons) invalid longitudes")
            else
                push!(messages, "✓ All longitudes valid")
            end

            # Compute statistics
            innovations = [obs.ob_minus_background for obs in info.observations]
            if !isempty(innovations)
                mean_innov = mean(innovations)
                rms_innov = sqrt(mean(innovations.^2))
                push!(messages, "✓ Mean innovation: $(@sprintf("%.4f", mean_innov))")
                push!(messages, "✓ RMS innovation: $(@sprintf("%.4f", rms_innov))")
            end
        end

    catch e
        push!(messages, "ERROR: Failed to read file: $e")
        is_valid = false
    end

    return (is_valid, messages)
end

"""
    compare_diagnostic_files(file1::String, file2::String; tolerance::Float64=1e-6)

Compare two diagnostic files and report differences.

Returns a comparison report as a vector of strings.
"""
function compare_diagnostic_files(file1::String, file2::String; tolerance::Float64=1e-6)

    println("\n" * "="^80)
    println("COMPARING DIAGNOSTIC FILES")
    println("="^80)
    println("File 1: $file1")
    println("File 2: $file2")
    println()

    report = String[]

    try
        info1 = read_diag_obs_file(file1)
        info2 = read_diag_obs_file(file2)

        # Compare file sizes
        size_diff = abs(info1.file_size - info2.file_size)
        size_pct = info1.file_size > 0 ? 100.0 * size_diff / info1.file_size : 0.0

        push!(report, "FILE SIZE COMPARISON:")
        push!(report, "  File 1: $(info1.file_size) bytes")
        push!(report, "  File 2: $(info2.file_size) bytes")
        push!(report, "  Difference: $size_diff bytes ($(@sprintf("%.2f", size_pct))%)")

        if size_pct > 10.0
            push!(report, "  ⚠ WARNING: File sizes differ by more than 10%")
        else
            push!(report, "  ✓ File sizes similar")
        end
        push!(report, "")

        # Compare record counts
        push!(report, "RECORD COUNT COMPARISON:")
        push!(report, "  File 1: $(info1.n_records) records")
        push!(report, "  File 2: $(info2.n_records) records")
        push!(report, "  Difference: $(abs(info1.n_records - info2.n_records))")

        if info1.n_records != info2.n_records
            push!(report, "  ⚠ WARNING: Record counts differ")
        else
            push!(report, "  ✓ Record counts match")
        end
        push!(report, "")

        # Compare observation counts
        push!(report, "OBSERVATION COUNT COMPARISON:")
        push!(report, "  File 1: $(info1.n_observations) observations")
        push!(report, "  File 2: $(info2.n_observations) observations")
        push!(report, "  Difference: $(abs(info1.n_observations - info2.n_observations))")

        if info1.n_observations != info2.n_observations
            push!(report, "  ⚠ WARNING: Observation counts differ")
        else
            push!(report, "  ✓ Observation counts match")
        end
        push!(report, "")

        # Compare statistics if both have observations
        if info1.n_observations > 0 && info2.n_observations > 0
            stats1 = compute_innovation_statistics(info1.observations)
            stats2 = compute_innovation_statistics(info2.observations)

            push!(report, "STATISTICS COMPARISON:")

            for key in ["mean_innovation", "rms_innovation", "bias", "std_innovation"]
                val1 = get(stats1, key, 0.0)
                val2 = get(stats2, key, 0.0)
                diff = abs(val1 - val2)

                push!(report, "  $key:")
                push!(report, "    File 1: $(@sprintf("%.6f", val1))")
                push!(report, "    File 2: $(@sprintf("%.6f", val2))")
                push!(report, "    Difference: $(@sprintf("%.6e", diff))")

                if diff > tolerance
                    push!(report, "    ⚠ WARNING: Values differ by more than tolerance")
                else
                    push!(report, "    ✓ Values match within tolerance")
                end
            end
        end

    catch e
        push!(report, "ERROR: Comparison failed: $e")
    end

    # Print report
    for line in report
        println(line)
    end

    println("="^80)

    return report
end

end # module DiagnosticReaders
